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Abstract. We study the properties of charmonium states at finite temperature in quenched 
lattice QCD on large and fine isotropic lattices. We perform a detailed analysis of charmo- 
nium correlation and spectral functions both below and above T c . Our analysis suggests 
that the S wave states disappear at about 1.5 Tc. The charm diffusion coefficient is esti- 
mated and found to be approximately 1/nT at 1.5r c < T < 3T C . 
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In the hot and deconfined QCD medium it has been conjectured that, unlike light mesons, the heavy 
mesons might survive and dissociate only at some higher temperature due to Debye screening (TJ. 
Therefore the suppression of their yield in nucleus-nucleus compared to proton-nucleus or proton- 
proton collisions may serve as probe for the properties of the medium. The experiments carried out 
at the SPS and LHC at CERN and the RHIC at BNL have indeed observed J/ifi suppression J2|. The 
interpretation of experimental data, however, is not as straightforward and in order to disentangle the 
cold nuclear matter effects on the one hand and hot medium effects on the other, it is crucial to have a 
good understanding of the behavior of heavy quarks and quarkonia in the hot medium. In addition, a 
substantial elliptic flow of heavy quarks has been observed [ 3 ] . As the heavy quark diffusion coefficient 
D can be related to the ratio of shear viscosity to entropy density r\j s I4I5II , a non-perturbative estimate 
of D is highly desirable. 

From the theoretical point of view, the meson spectral function at finite temperature contains 
all the information on the hadron properties in the thermal medium. This explicitly includes the infor- 
mation on bound states and transport properties. As a consequence the non-perturbative computation 
of spectral functions is of key interest. 

In the following we report on a recent success in this direction in quenched lattice QCD. We present 
the spectral functions in the pseudoscalar and vector channels and estimate the heavy quark diffusion 
constant. See [7| for a recent overview and [8| for a more detailed and extended discussion. 
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Table 1. Lattice parameters and number of configurations used in the analysis with a clover improved Wilson 
fermion action. 





a[fm] a"'[GeV] 


U [fm] 


csw 


K N^XNr 


T/Tc 


conf 


7.793 


0.010 18.974 


1.33 


1.310381 


0.13200 1 28 3 x 96 


0.73 


234 










128 3 x48 


1.46 


461 










128 3 x 32 


2.20 


105 










128 3 x 24 


2.93 


81 



2 Euclidean correlation functions and the charmonium spectral function 

At vanishing momentum charmonium correlation functions are defined via current operators by: 

Gh(t, T) = J]< JhM 4(0,0) ) T . (1) 

X 

Here Jh is e.g. a local mesonic current operator q(j, x)rH<l(T, x), where = yujs for vector (V„) 
and pseudo-scalar (PS) channels, respectively. Notice (...)r denotes a thermal average. On the lattice 
the temperature T is related to Euclidean temporal extent aN T by T — 1 / (aN T ), where a is the lattice 
spacing. The above correlation function is exactly related to the spectral function via: 

r*oo J*. 

G h (t,T)= — p H (o,T)K(T,T,aS), (2) 
Jo 2tt 

where the kernel K is given by A^t, T, to) = cosh(<w(T - ^))/sinh(^). The spectral function p(u>) 
contains all the information on the in medium properties of hadrons and its calculation therefore is 
our primary goal. The dissociation temperature can be read off from the deformation of the spectral 
function, while the heavy quark diffusion constant D is related to the vector spectral function via the 
Kubo formula: 

D = J- (3) 

6*00 (*> 

where *oo is the quark number susceptibility. 

To study and estimate the charmonium spectral functions in the pseudoscalar and vector channels 
as temperature is increased, we first examine the behavior of the correlator itself in a number of ratios 
and differences. It should be clear that in this fashion we can only obtain a qualitative idea on the 
properties of the spectral function. The estimates are then used to construct the necessary default 
models for an analysis based on the maximum entropy method (MEM) (5). Using MEM we invert 
Eq.[2]and compute the desired spectral function. A detailed analysis of our systematics and statistical 
uncertainties can be found in |8] . 



3 Numerical results 



For our numerical investigation the standard Wilson plaquette action was used in the gauge sector, 
while for the (valence) charm quarks we choose the non-perturbatively 0(a) improved Wilson clover 
fermion action. Note, here the mass of vector meson is tuned to be almost the physical J /iff mass. 
Employing a fixed scale approach we measure correlation functions on very fine quenched lattices 
sized 128 3 x 96, 128 3 x 48, 128 3 x 32 and 128 3 x 24. These correspond to the temperatures 0.73 T c , 
1.46 T c , 2.20 T c and 2.93 T c , whereby T c * 270MeV in quenched QCD. Note, the required lattice 
parameters are given in Tab. Q] 
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Fig. 1. Euclidean correlation functions G(r)[GeV 3 ] in the pseudoscalar PS (left) and vector V,,(right) channels at 
P = 7.793 (a = 0.010 fm) and the temperatures T = 0.73, 1.46, 2.20 and 2.93 T c 



3.1 Charmonium correlation functions 

To begin our analysis we show the correlation functions in the pseudoscalar and vector channels in 
physical units in Fig. [TJ Clearly all correlators agree throughout most of the available Euclidean time 
distance except when they near their respective midpoints. This small effect in the large distance re- 
gion is the essential piece for understanding the thermal modification of the bound state and transport 
regions of the spectral function, while the correlation function at small to intermediate distances is 
dominated by the continuum contribution. As a consequence a more sophisticated analysis highlight- 
ing the thermal modifications by taking e.g. appropriate ratios or differences is necessary. 

To analyze the properties of the underlying spectral functions on the correlator level more thor- 
oughly we turn to studying it in conjuncture with the so called reconstructed correlator. 
In a reconstructed correlator the temperature dependence of the kernel in Eq.|2]is taken into account, 
while the spectral function remains unchanged. It is particularly interesting to compare a zero or low 
temperature spectral function with a thermal one and as it turns out in this case it is possible to calculate 
the reconstructed correlator analytically as follows, first one constructs 1101 : 

f°° dco , cosh(w(r- 1/2T)) 

G rec (r,T;T') = — p(u, T') \ V . ' (4) 

Jo 2n sinh(t<j/2r) 

At this point one would generally need the spectral function p(a>, T') at a reference temperature 
T and consequently the evaluation of Eq. (|4]i would suffer from the uncertainty of the determination 
of that spectral function. However there is a useful and exact relation, which is a generalization of 
Ref. IfTTI : 

coshMr - N T /2)] K y +T coshMV ~ K/2)] 

sinh(wjV T /2) = 2-j sinh(wjV;/2) ' 

where V = (aAQ -1 , T = (aiV T ) _1 , r' e [0, N' T - 1], r e [0, N T - 1], N' T = m N T , m e Z + . N T 
and N' T are the number of time slices in the temporal directions at temperature T and T', respectively. 
t denotes the time slice of the correlation function at temperature T while t' denotes the time slice of 
the correlation function at temperature T . The sum of r' on the right hand side of Eq. (0 starts from 
t' — t with a step length A^ T to the upper limit N' T — N T + t. After putting p(co, T') into both sides of the 
above relation and performing the integration over a>, one immediately arrives at If8l : 

N' t -N t +t 

G rec (T,T;T') = J] G(r',r'), (6) 

t'=t; At'=N t 

with this relation we are therefore capable of computing the reconstructed correlator directly from our 
low temperature results at A^ T = 96 without having to extract the spectral function, thereby bypassing 
a systematic source of error in this quantity. 
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Fig. 2. The ratio G(t, T}/G rec (r, T) for PS (left) and V„ (right) channels as a function of the Euclidean distance 
t at J = 1.46, 2.20 and 2.93 r c . The reconstructed correlator G rec is obtained directly from correlator data at 
0.73 T c . 



Given this insight we use Eq.|6]to compute the ratio of the thermal and reconstructed correlators 
G(j)IG rec {T) at 1.46r c , 2.207\. and 2.93T C in the pseudoscalar and vector channels, where the results 
are shown in Fig. [2] Concentrating first on the pseudoscalar case in Fig. |2|left) we observe that the 
ratios decrease monotonically at the available temperatures as the Euclidean time separation increases. 
Note, at the largest distances, the ratios deviate from unity by about 5%, 8% and 12% at 1.46, 2.20 
and 2.93 T c , respectively. Turning to the right of Fig.|2]and therefore the V„ channel we see an entirely 
different behavior. Here the ratios increase monotonically with increasing distance, reaching a devia- 
tion of about 16% at the large distances. These results are interesting as there is no emerging transport 
contribution in the PS channel, therefore our results hint at the decrease in this channel being due 
to bound state modification, while a large part of the modification in the V„ channels is due to the 
emerging transport contribution. 

To study this more closely we define two ratios, which are approximately independent of a zero 
mode contribution: 

G diff (r, T) = G(t,T)-G(t+1,T) G sub (T, T) = G(t, T) - G(N T /2, T) 

G^(t, T) = G rec (T,r)-G rec (T+l,r)' an G^(t, T) = G kc (t,T)-G kc (N t /2,TY 

The first equals the ratio of the time derivative of the correlators to the time derivative of the recon- 
structed correlators at t+ 1/2 [ 12|, while the second is the ratio of midpoint subtracted correlators [ 13 1. 
As noted above both ratios are approximately independent of a zero mode contribution, consequently 
they should be almost independent of transport phenomena. The corresponding results in the PS chan- 
nel are given in Fig.[3]left) and we find that both G fl?, ^(r)/G*{^(r) and G™ i> (T)/G™*(r) give very sim- 
ilar results. Throughout the available temperature range we retain the decreasing trend with increasing 
distance from the analysis of G(T)/G rec (T), although the magnitude of the deviation has decreased, the 
values at the largest distance being shifted up by about 3% at both 1.46 T c and 2.20 T c and about 6% 
at 2.93 T c . As there is no transport contribution in this channel we would have expected fairly small 
changes from the behavior of G(T)/G rer (r) to G d 'f J \t) I G d r ' e { 3 \t) or G ml '(t) / 'G™^(t) and indeed we see 
this confirmed. 

Turning to the V u channel on the right of Fig.|3]we first note that G dif f (t) I Gf/J (t) and G sub (f)/G s r £(T) 
give similar results as in the PS case. The behavior compared to G(r) / G rec (T) on the other hand has 
been greatly modified. As such the data for 1.46r e instead of increasing now is almost flat, while 
2.20r e and 2.93r c show a decreasing trend similar to the PS case. From this we infer that most of the 
temperature dependence in G(T)/G ree (r) in the Vg channel is in fact due to the transport contribution. 
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Fig. 3. The ratio G dlff (r, T)/G;™(t, T) (G sllb (T, T)/G™ b (T, T)) for the PS (left) and V,,(right) channels as a function 
of the Euclidean distance t at T = 1.46, 2.20 and 2.93 T c . The superscript "x" denotes either "diff' or "sub". 



Another way to look at the data is to take the difference G(t) - G rec (T), implying: 

— 

o 2tt 

consequently this difference between the measured correlator and the reconstructed correlator provides 
information on the difference of spectral functions below and above T c and the corresponding results 
are given in Fig. |4] Note in both V,,- and PS channels we observe that G(r) - G rec {T) increases with 
tT at the available temperatures. In the PS case on the left of Fig. [4] however all values are negative, 
while for the V„ on the right they are positive. 

As G(r) - G rec (t) probes the difference of the two spectral functions Ap(uS) the results in the PS chan- 
nel require Ap{io) < in a certain to region. With no transport contribution in this channel and given 
the positivity condition of the spectral functions this implies the low temperature bound state peak 
overcompensates any remnants of the peak at finite temperature. 

At the same time in the V„ channel the positive values indicate a more intricate situation, as the emer- 
gence of transport leads to a positive contribution, while an overcompensation of the bound state would 
lead to a negative part. 

Note, the weakest r-dependence is seen in the V„ channel at 1 A6T C . Here we observe an almost flat be- 
havior with positive values throughout. Assuming that there is a complete cancellation of bound state 
peaks at this temperature in Ap{to) we may try to fit the transport contribution to the Ansatz defined in 
Eq. [3] In this case there is only one free parameter and the uncertainty in the heavy quark mass, for 
M = l.OGeV and M = 1.8GeV we obtain: 

M = l.OGeV, 2kTD « 0.6 and M = 1.8 GeV, 2nTD * 3.6. (9) 

Even though it is not convincing to use a one parameter fit here to obtain a robust value for the heavy 
quark diffusion coefficient, these results are in the right ballpark and provide useful input for our 
subsequent MEM analysis. 



p(a>,T)-p(a), 7") 



cosn(Mr - i/^i )) 
sinh(w/2D 



(8) 



3.2 Charmonium spectral functions via MEM 

Clearly the comparison of the measured correlator with the reconstructed correlator can only give a 
rough idea of the magnitude of any medium effects at a certain temperature. To really explore the 
properties of the thermal charmonium states, one has to go to the level of spectral functions. 
Here we do this by using the maximum entropy method (MEM) to extract the spectral functions 
directly from the correlator data. In the MEM analysis presented here, we use as number of fre- 
quencies Naj = 8000, the minimum energy aa> mi „ = 0.000001 and implement the modified kernel 
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Fig. 4. (G(t) - G,-„ (t))/T 3 in the .PS (left) and V,,(right) channels as a function of the Euclidean distance t at 
T = 1.46, 2.20 and 2.93 T c . 




Fig. 5. Statistical uncertainties of output spectral functions in PS (left) and V„ (right) channels at all available 
temperatures. The shaded areas are statistical errors of amplitudes of output spectral functions from Jackknife 
analyses and the solid lines inside the shaded areas are mean values of spectral functions. The horizontal error 
bars at the first peaks of spectral functions at 0.73 T c an 1.46 T c stand for the statistical uncertainties of the peak 
location obtained from Jackknife analyses. 



K* = tanh(d)/2) • K H 141 1 51 in order to explore the low energy behavior of spectral function [ 16 1. The 
default model in the PS channel is a normalized free lattice spectral function, while in the V„ channel 
it is a normalized free lattice spectral function with a transport peak based on Eq.[3]and the diffusion 
coefficient estimated above. 

The results of this analysis are given in Fig. |5J where we show the PS (left) and y,,(right) spectral 
functions divided by to 2 at the available temperatures above T c and also below T c . Note, through- 
out we show the spectral functions in PS and V„ channels with error bands given by the statistical 
uncertainties, while the error of the peak positions, see (SJ, are given as horizontal error bars. 

In Fig.|3left) we observe that at 0.73 T c the spectral function in the PS channel has large uncer- 
tainties in the amplitude at the point which corresponds to the ground state peak location in the mean 
spectral function. However, even at the lower end of the error band, the amplitude is still larger than 
the peak amplitudes at the higher temperatures within the errors. 

Unlike the large uncertainties shown in the amplitude of the peak height, the peak location of the 
ground state peak at 0.73 T c is well determined. A Jackknife analysis yields m 1tc = 3.31(4) GeV. At 
1.46 T c this peak is shifted by about 0.8 GeV to around 4.1 GeV. At 2.23 T c there is hardly a peak 
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Fig. 6. Left: statistical uncertainties of transport peaks at T > T c . Right: the resulting charm diffusion coefficients. 
The boxes stand for statistical error estimated from Jackknife method while the bars stand for systematic uncer- 
tainties from MEM analyses. To compare also the results from perturbation theory 1171 . Langevin dynamics (3), 
Ads/CFT COD and HQET I2TT221 are given. 



structure that can be identified within the statistical uncertainties. At 2.93 T c the spectral function 
flattens further. Thus this picture suggests that j] c is melted already at 1.46 T c , 

Turning to the right of Fig. [5J we analyze the resonance part of the spectral function in the V„ 
channel. We observe that the peak location of the spectral function at 0.73 T c does not have an overlap 
with the peak location of the spectral function at 1 .46 T c and the amplitudes between these two differ a 
lot (see horizontal error bars). At both 2.20 T c and 2.93 T c there is hardly any peak structure. As such 
this picture indicates that also J/ift is already dissociated at 1.46 T c . 

Focusing on the very low frequency part of the spectral function given in the vector channel, we 
show the spectral function divided by coT in Fig. |6lleft). Notice here the statistical uncertainties on the 
amplitude of the peak are relatively small. The charm diffusion coefficient is related to the amplitude 
of the transport peak at vanishing frequency through the Kubo formula (Eq. [3). The corresponding 
estimates for the heavy diffusion coefficient D are summarized in the right plot of Fig. [6] The boxes 
denote the statistical uncertainties and the error bars reflect systematic uncertainties obtained from the 
analyses discussed in [8 |. Note the charm diffusion coefficient obtained at 1.46 T c is the most reliable 
one among the three temperatures above T c since more prior information is known at this temperature. 
At 2.20 T c and 2.93 T c , due to the lack of precise prior information and also less number of data points 
that can be used in the MEM analyses, the uncertainties on the charm diffusion coefficient thus might 
be underestimated. 

For comparison we highlight results from alternative approaches to the calculation of the heavy 
quark diffusion coefficient in Fig. fright), showing the perturbative next-to-leading order results of 
[17 1, those from a computation in AdS/CFT [18], as well as those from Langevin dynamics derived 
in 0. Additionally we show results from a recently proposed method to compute the momentum 
diffusion coefficient k from a purely gluonic operator derived from heavy quark effective theory 
(HQET) [19 20 1 . Hereby k is related to the heavy quark diffusion coefficient by k — 2T 2 /D. First 
calculations from two groups on this operator agree fairly well [21 .22] and compare favourably with 
our results. 



4 Summary 

In summary our analysis of charmonium correlation functions on large isotropic quenched gauge con- 
figurations using Wilson clover fermions implies that both J/ifr and rj c have dissociated by a tempera- 
ture of 1.467V At the same time we estimated the charmonium diffusion coefficient and found a value 
2nTD ^ 1.8 - 2.3 in the available temperature range. 
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